Vortex Core Structure in Neutral Fermion Superfluids with Population Imbalance 
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Quantized vortex core structure is theoretically investigated in Fermion superfluids with popula- 
tion imbalance for two atom species of neutral atom clouds near a Feshbach resonance. In contrast 
with vortex core in balance case where the quantum depletion makes a vortex visible through the 
density profile measurement, the vortex core is filled in and becomes less visible because the quan- 
tized discrete bound states are occupied exclusively by the majority species. Yet it is shown that the 
core can be visible through the minority density profile experiment using phase contrast imaging, 
revealing an interesting opportunity to examine low-lying Fermionic core bound states unexplored 
so far. 

PACS numbers: 03.75.Ss,03.75.Hh,47.32.-y 



There has been much attention focused on Fermionic 
superfluids of cold atoms, such as 6 Li or 40 K since ex- 
perimental realization The experiments utilize the 
Feshbach resonance by changing an external field to 
tune the atom-atom interaction, achieving the BEC-BCS 
crossover. A keen interest is now placed on Fermion su- 
perfluidity when the population of the two species (up 
and down "spins") is unequal experimentally [3, d, 0, Q 
and theoretically p, 0] ■ Zwierlein et al. 0, Q have suc- 
ceeded in creating vortices not only in the balance case 
with equal population, but also in the imbalance case, di- 
rectly demonstrating its superfluidity. Here the presence 
or absence of vortices in a system is utilized to monitor 
superfluidity because quantized vortex is a hallmark of 
superfluidity. They have observed the clear signature of 
vortices in the inner region of the density after sweeping 
into the BEC side, even in the imbalance case of the BCS 
side. However, it is important to notice that whether the 
superfluidity in the outer region is robust or not is still 
open to question. In this Letter, we focus on the visibil- 
ity of the vortex core situated at the trap center in the 
imbalance case. In addition, we discuss the visibility of 
vortices in the outer region from the results. This result 
is quite contrasted in the BEC case where vortices are 
remarkably arranged regularly throughout whole system, 
even near the boundary with lower density 0]. 



The order parameter is ^bec = \/n(r) (n(r) is the 
number density) in BEC, while the Cooper pair ampli- 
tude ^bcs oc (c£jcl_ fc |) in BCS. Thus it is not self-evident 
even in the 50%-50% balance case that the vortex, which 
is probed by the density contrast, continues to be "seen" 
from BEC to BCS across the Feshbach resonance point. 
A theoretical question here is to understand the peculiar 
vortex core structure in the imbalance case, related to 
the quantum depletion of the density at the core 0, [l(| • 
We provide a microscopic calculation for it, fully taking 
account of low-lying Fermionic excitations around a core 
beyond simple local density approximation. These indi- 
vidual excitations are essential in Fermion superfluid, a 
feature completely absent in the vortex state in Boson 



superfluid. 

This study might be useful for other research fields 
such as condensed matter community because super- 
conductors with imbalance up and down spin electrons 
are realized by an applied field through the Zeeman ef- 
fect. Thus there is a good chance to check the present 
study. Indeed in a heavy Fermion material CeCoIn 5 the 
"imbalance" superconductivity may be realized under 
high fields where the Fulde-Ferrell-Larkin-Ovchinnikov 
(FFLO) state is observed [HI]. Also in high temperature 
superconductors the so-called "vortex charging" problem 
where the quasi-particle density is suppressed at the core 
as observed by nuclear quadrupole resonance (NQR) ex- 
periment 12 1 and anomalous Hall conductivity [13[. We 
also expect that the study of the vortex core in the imbal- 
ance superfluid might be useful for the inside of rotating 
neutron stars where vortices are sustained in imbalance 
color superconductivity, consisting of quark matter (l3 |. 

Before going into detailed calculation in the imbalance 
case, we recapture the results in the balance case; The 
core structure both in superconductors and superfluids 
is studied in the balance case [t| [l(| • They show that as 
increasing the coupling constant of attractive interaction, 
the particle number at the core decreases due to quantum 
depletion. This is caused by the fact that the vortex 
core bound states, namely, Caroli-de Gennes-Matricon 
(CdGM) states [l|| with finite amplitude at the vortex 
center are unoccupied by discretization |T(| • It should 
be noted, however, that there is no experiment to directly 
detect the CdGM state in condensed matter physics so 
far. Here we show through this study that there is a good 
chance to directly see it. 

To address such a problem, we self-consistently solve 
the Bogoliubov-de Gennes (BdG) equation describing the 
interacting Fermion system with population imbalance 
under a trap V(r). This framework is one of the most 
fundamental microscopic theories. We start with the 
BdG equation for the quasi-particle wave functions u q (r) 
and Vq(r) labeled by the quantum number q as follows 
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where K u (r) = + V"(r) + sni, T (r) - Mu- Wc 

fully takes account of the mismatch Fermi surfaces as 
Sfi = (/if — /Ui)/2 and the local Hartree potential 
gn a {r) with n T (r) = £ q |u q (r)| 2 /(£q) and n i( r ) = 
S q l^q( r )| 2 [l — /C^q)]- While eigenstates and eigen- 
value of Eq. (1) are often written as (u q ',|, Wq',t)> 



Eqr^, we express these states 



-Eq'.l and (- 

as (u q ,Wq), i? q . The bare attractive interaction g 
Awa/M via the scattering length a. Throughout this 
paper, we set h = 1 and consider a cylindrical sys- 
tem with V(r) = |MwV, imposing a periodic bound- 
ary condition with the periodicity Z = 3d {d = 
yMcij) along the z-direction. Hence, we write the 
wave functions as u q (r) = u q (r) exp [i(£ — i)6 + iq z z] 
and Vq(r) = u q (r) exp + + iq z z] with q z = 
0,±27r/Z,±47r/Z, ■ • • for a single-vortex. The BdG ma- 
trix in Eq. (1) is numerically solved by spatial discretiza- 
tion with respect to r-axis, whose step is O.Old. The 
self-consistent gap equation for A(r) = e~ l9 A(r) is 



A(r) =g(r) 



E 

|E q |<B c 



WqM^qM/CEq), 



(2) 



where the effective coupling constant g(r) is introduced 
in order to remove the ultra-violet divergence [H, 
Since we use the expression valid even for broken particle- 
hole symmetry, the sum in Eq. (2) is done for all the 
eigenstates with both positive and negative eigenenergies 
up to the energy cutoff ±E C = 200w [11]. In the BdG 
equation (TTJ), the chemical potential is adjusted to fix 
the total particle number TV = 7V T + 7V X = 3000. Here, 
we shall present the self-consistent results for the bare 
coupling constant g/wd 3 = —1.1, —1.2, —1.3, —1.4, —1.5, 
corresponding to —1.4 < (kpa)^ 1 < —0.7 with the 
Fermi wave number hp. These are also expressed as 
A (0)/£ F (0) = 0.15,0.20,0.23,0.27,0.32 in the vortex- 
free state, respectively. The gap A o (0) (Ep(Q)) is A(r) 
(the Fermi energy Ep) at the trap center without a vor- 
tex in the population balance case. The polarization P 
is defined by P = (iV T - iV;)/(iV T + JV;). In this paper 
we only consider the T — case and neglect the external 
rotation energy. 

We show the spatial profiles of the order parameter 
A(r) in Fig. 1. The vortex structure for the population 
imbalance corresponds to curve (1). It is seen that the 
vortex core size is widen compared with that for the bal- 
ance case (2). Towards the outer region A(r) oscillates 
and changes its sign several times to accommodate the 
excess majority component. This FFLO-like state is al- 
ways stabilized irrespective of the values of finite P at 
T = as shown previously @. Therefore for P ^ 



at T = in the vortex solution A(r) changes its sign 
and oscillates. As shown in the inset, which displays 
|A(r)| for curve (1) in a logarithmic scale to examine 
the detailed variation of the outer region, the oscillation 
pattern continues rather regularly to smaller scales of A. 
This means that superfluidity is kept coherently through- 
out the whole system, even at the system edge where |A| 
is extremely small and the minority component is almost 
absent (See Fig. 2). Therefore the present solution for 
the imbalance case is quite different from a simple BCS- 
normal phase-separation in the local density approxima- 
tion. This is consistent with that in our previous solution 
for the vortex free case Q and with the work by Kinnunen 
et al. 0]. 




FIG. 1: (Color online) Order parameter profiles |A(r)| as a 
function of r/d for g/iud 3 = —1.5. (1) Vortex with population 
imbalance case P = 0.3. (2) Vortex with population balance 
case P = 0. (3) Vortex free with population balance case 
P = 0. Logarithmic plot of (1) is shown in inset to emphasize 
the finer scale of |A(r)| in the outer region. 

In Fig. 2 we show the density profiles; the total n(r), 
the majority n-\(r) and minority n^(r) components, and 
the magnetization m(r) — Jij(r) — n\[r) for the solution 
corresponding to curve (1) in Fig. 1. It is seen that the 
majority component exclusively occupies the central re- 
gion while the minority component remains depleted. It 
should be noted that the quantum depletion of the total 
number n(r) at T = is substantially smaller and in- 
complete compared with the balance case as seen shortly 
in Fig. 4. The magnetization shows a peak at the vortex 
center. In the intermediate region 0.5 < r/d < 2.5 the 
magnetization is almost vanishing, resulting in the en- 
hancement of ni(r). This leads to the "bimodai" struc- 
ture of ni(r), which is observed experimentally in the 
vortex free cases 0, H[- This bimodai structure is con- 
firmed by our previous calculation in the vortex free case 
Q . In the outer region r > 5d in Fig. 2 there exists only 
the majority component n-\ (r) and the magnetization is 
fully polarized there. 

In Fig. 3 we show density plots for n^{r),ni(r) and 
n(r) of Fig. 2 where the vortex core is located at the 
center. The contrast in ny{r) and n(r) is clearer than 
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FIG. 2: (Color online) Density and polarization profiles as a 
function of r/d corresponding to curve (1) (P = 0.3, g/ojd 3 = 
— 1.5) in Fig. 1. n(r) — rif(r) + rij(r): total density, 
nt(r)(ni(r)): up (down) spin density and m(r): magneti- 
zation. 

that of n-[(r), which is dim. Namely in the imbalance 
vortex, the minority (majority) component is visible as 
in Fig. 3(b) (invisible as in Fig. 3(a)). We also show the 
density plot of n(r) for the balance case for comparison in 
Fig. 3(d). It is seen that the contrast of the balance case 
is clearer compared with that for the imbalance case in 
Fig. 3(c). This is because the quantum depletion of n(r) 
occurs at the core both for n^(r) and ri|(r). These results 
in Fig. 3 are checked experimentally by direct "in situ" 
method [H, 0] without resorting to BEC mapping from 
the BCS side, combined with phase contrast imaging 
When the BEC mapping is used [2| , the density contrast 
of vortices may be modified. We note in the passing that 
the BEC case n(r) = at the vortex center, because n(r) 
is the order parameter and at the vortex center n(0) must 
vanish because of the presence of the phase singularity 
due to phase winding. 

The polarization P dependence of the total density 
profiles is displayed in Fig. 4. In the balance case (P = 0) 
substantial quantum depletion at the vortex core is seen. 
As P increases, the core density is progressively filled out 
and the vortex filling becomes complete, thus the vortex 
tends to be "invisible" in the density contrast experi- 
ment. The core radius defined at the density maximum 
point from the trap center increases with increasing P, 
then the "core size" expands. 

The inset of Fig. 4 shows the P-dependence of the core 
filling factor F = n(0)/n max defined by the density ra- 
tio between the density n(0) at the core center and the 
maximum density n max . It demonstrates how visibility 
becomes lower with P. Due to the quantum depletion 
at the core, which is proportional to Ao(0)/Pf [S|, the 
filling factor F becomes smaller as the strength of the 
attractive interaction \g\ increases, i.e. approaching the 
Feshbach resonance point. The density contrast quickly 
diminishes with P, implying that vortex becomes invis- 
ible with P for general g. As seen from the inset of 




FIG. 3: Density maps of (a) n-\(r), (b) n;(r) and (c) n(r), 
corresponding to Fig. 2. (d) n(r) for P = and g/u)d 3 — —1.5. 
By comparing the imbalance case (c) and the balance case (d), 
the core in (d) is clearly seen as a small but sharp dot at the 
center. The field of view is 14.0 x 14.0 in units of d. 
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FIG. 4: (Color online) Density profiles of the total number 
n(r) for various polarizations: P = 0,0.11,0.35,0.43,0.52 
from top to bottom. The inset shows the filling factor 
F = n(0)/n max defined by the ratio between the maximum 
value n m ax and n(r = 0) for various couplings {g/tod 3 = 
-1.5, -1.4, -1.3, -1.2, -1.1 from bottom to top). 

Fig. 4 the filling factor F as a function of P changes 
in a step-like manner, particularly for larger \g\. This 
comes from the fact that the vortex bound states are 
discretized. Thus as P varies, the occupation of these 
quantized states is progressively changed, leading to the 
step-like variation of the occupation number. 

These results might be relevant for the interpretation 
of the vortex experiment in the imbalance case, espe- 
cially, around the outer region. In the actual experi- 
ment the vortex visibility is influenced by the formation 
of FFLO as follows: (i) In the outer region, FFLO pairing 
survives with small amplitude and (ii) it is also always 
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accompanied with a large amount of polarization. By 
these effects, vortices around outer region in the FFLO 
state can be invisible. The vortices at the center (outer) 
correspond to the present single vortex study for smaller 
(larger) P. Hence, we expect general tendency that in the 
imbalance case the vortices in the center (outer) region 
are visible (invisible). 
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FIG. 5: (Color online) Spectral evolutions for the majority 
species (up-spin atoms) Nf(r,E) with q z = 0. Stereographic 
view (a) and density map (b). P — 0.3 and g/ud? = —1.5. 
The labels show the angular momentum I for each CdGM 
bound state. In (b) the order parameter ±A(r) is overlaid. 
Energy is normalized by the harmonic frequency u. Ef is the 
Fermi level position. 

In order to understand why the total number n(r), 
in particular the majority species n^(r) fills out the 
core region, we examine the local density of states with 
q z = 0, namely, JV T (r,£?) - £ q \ Utl (r)\ 2 S(E - £ q ) and 
N l( r , E ) = E q \vq{r)\ 2 6(E+Eq). The spectral evolution 
of Nf(r,E) for the majority species near the Fermi level 
Ep and the corresponding density plot are displayed re- 
spectively in Figs. 5(a) and 5(b). It is seen clearly from 
Fig. 5(a) that near the vortex core r = and at the 
lower energy region around E — 0, a series of the core 
CdGM bound states 15j, i.e. more broadly the Andreev 
bound states, are located, whose angular momentum is 
characterized by half- integers £ = ±i,±f,---. At the 
core center r = only the CdGM bound state labeled 
by t — I with a positive energy has spectral amplitude 
and the others have no spectral weight at r = [9] . This 
bound state is unoccupied when P = 0. Since for P ^ 



it is occupied because of the upward Fermi level shift, 
the majority species n^(r = 0) becomes able to fill in the 
core. For the minority species this state is unoccupied, 
thus n j (r) remains suppressed at the core. Since the high 
angular momentum states with I = 3/2, 5/2, • • • are pro- 
gressively occupied for the majority species, n^(r) further 
increases near the core. Beyond a certain P value, which 
depends on the coupling constant, the depletion of the 
core becomes absent and the core filling is complete (see 
inset of Fig. 4). 

We should notice from Fig. 5(a) that each CdGM 
bound state is split and has the double peak structure. 
This can be understood as follows: We show the spectral 
density map of N* (r, E) together with the order param- 
eter profile ±A(r) in Fig. 5(b). At the node (r/d = 3.0) 
where A(r) = 0, there appear the zero-energy bound 
states, i.e. the so-called midgap state associated with 
FFLO. This midgap state accompanies several excited 
states around the node because quasi-particles feel a lo- 
cal confining potential created by the nodes of |A(r)| 
[2fj| . Although these Andreev bound states are local- 
ized around the nodes, they can resonate with the CdGM 
state with £ = | bounded by the vortex if these ener- 
gies match each other. Other CdGM states with higher 
angular momentum may also resonate with the surface 
harmonic states which are originally localized at the sur- 
face. These new "split" bound states further enhance the 
n-f occupation and the magnetization around the vortex 
core. Physically this enhancement is reasonable because 
these newly created resonant bound states are advanta- 
geous in view of the potential energy gain at the center. 

In summary, we have shown by solving microscopic 
BdG equation self-consistently that the vortex core in 
Fermion superfluids tends to be filled in by the majority 
species when the population of the two species is im- 
balance. It makes the vortex core image unclear via the 
density profile measurement. The vortex core can still be 
visible through the selective measurement of the minor- 
ity density, using phase contrast imaging. We have given 
its physical reasons in terms of the Caroli-de Gennes- 
Matricon bound states. These predictions can be checked 
experimentally by direct "in situ" method [!, 0| without 
resorting to BEC mapping Q from the BCS side, com- 
bined with newly developed phase contrast imaging 
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